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We perform large scale finite-temperature Monte Carlo simulations of the classical Cg and t2g orbital models 
on the simple cubic lattice in three dimensions. The eg model displays a continuous phase transition to an 
orbitally ordered phase. While the correlation length exponent v ~ 0.66(1) is close to the 3D XY value, the 
exponent -q ~ 0.15(1) differs substantially from 0(N) values. At Tc a (7(1) symmetry emerges, which persists 
for T < Tc below a crossover length scaling as A ~ with an unusually small a ~ 1.3. Finally, for the t2g 
model we find a. first order transition into a low-temperature lattice-nematic phase without orbital order. 

PACS numbers: ()5.70.Fh, 64.60.-i, 75.10.Hk, 75.40.Mg 



Orbital -only models emerged recently as prototype systems 
enabling the understanding of relevant aspects of the collec- 
tive dynamics of orbital degrees of freedom [ 1 ] . In a different 
context, orbital-like models are attracting considerable theo- 
retical interest due to their ability to sustain topologically or- 
dered phases with possibly anyonic excitations, as exemplified 
by the Kitaev honeycomb model [2] . In a similar spirit the or- 
bital compass model [3] can serve as a basic model to under- 
stand topologically protected Josephson junction qubits [4], 
which have recently been realized experimentally [5]. 

A variety of properties have already been uncovered for 
orbital-only models, but most of these are restricted to ground 
state or low-temperature properties. Much less is known about 
finite-temperature properties and in particular the nature of 
thermal phase transitions. Those might display new criti- 
cal phenomena, as a common feature of all these systems is 
a manifest coupling between order parameter space and real 
space, which distinguishes them from the well studied 0(N) 
(such as Ising, XY and Heisenberg) models [6]. 

In this Letter, we present a comprehensive Monte Carlo 
(MC) investigation of the nature of the finite-temperature 
phase transitions in two popular orbital-only models on the 
three-dimensional (3D) cubic lattice: the Cg and the t2g mod- 
els [1]. We study here the classical versions because the cor- 
responding quantum models have a sign problem precluding 
Quantum Monte Carlo approaches, and because in Ginzburg- 
Landau theory one typically expects quantum and classical 
versions of a same model to have the same critical proper- 
ties, although exceptions are possible. The eg and the t2g 
models are also often called the 120° and compass models, 
respectively. While the thermal phase transition in the two- 
dimensional (2D) compass model has been the focus of recent 
studies [7-9] and clarified to belong to the 2D Ising universal- 
ity class, little is known about the Cg and t2g models in 3D - al- 
though potentially of more direct relevance for the description 
of collective orbital phenomena [1]. We start by discussing 
the eg model and its critical properties in some detail and turn 
then briefly to the t2g model towards the end of this paper 



The eg model — The eg model (EgM) is defined by the 
Hamiltonian [1] 

^e, (1) 

where ti is an auxiliary three component vector obtained 
by an embedding of the orbital degree of freedom — 

/-1/2 V3/2 \ 
= -1/2 -73/2 T,. (2) 

V 1 y 

The Ba denote the positive unit vectors in the a £ {x, y, z} 
cartesian directions. Note that the coupling in x-space de- 
pends on the spatial orientation of the bond. The coupling 
constant J is set to one in the following, coiTesponding to 
ferromagnetic interactions. Note that results for antiferromag- 
netic interactions can be deduced from results using ferromag- 
netic couplings [10]. 

The classical EgM (1) has a sub-extensive ground state de- 
generacy which is lifted at finite temperature by an order by 
disorder mechanism [11], leading to six discrete ordering di- 
rections T° = (cos[n 27r/6], sin[n 27r/6]) with n = 0, . . . , 5. 
This analytical prediction has been verified using classical 
MC simulations [10, 12], and at higher temperatures a con- 
tinuous phase transition to a disordered phase has been found. 
The prominent question of the universality class of the finite- 
temperature phase transition is however still open, both ana- 
lytically and numerically. Comparing to related systems with 
a similar low-temperature phase, several different scenarios 
seem possible: i) a continuous transition in the universality 
class of the 3D XY model, as e.g. in the ^'g-perturbed XY 
models [13-16], ii) distinct universality classes, as reported 
in classical dimer models on the cubic lattice [17, 18] or iii) 
a first order transition, as in a six-state feiTomagnetic Potts 
model in 3D, or a Heisenberg feiTomagnet with a specific cu- 
bic anisotropy [6]. In the following, we shall resolve this fun- 
damental question and answer which scenario is realized for 
the eg and t2g models. 
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FIG. 1. (Color online), eg model: The order parameter m (left plot) 
and the associated Binder cumulant Bm (right plot) as a function of 
temperature T for different linear system sizes L. The vertical line 
indicates the location of the critical temperature. 



Simulation technique and observables — We consider the 
classical Hamiltonian (1) on a simple cubic lattice of side 
length L and volume N ^ and perform state-of-the-art 
MC simulations along the lines of Refs. [8, 9]. Simulations 
were performed for lattice sizes L = 8, . . . , 96. To obtain 
the reported accuracy, we collected 10^ and more independent 
MC measurements per data point. MC runs using periodic 
boundary conditions (PBC) show clear signals of a transition 
to an ordered phase in accordance with Ref. [10]. However, 
as further demonstrated below, we find that there are severe 
finite-size corrections using PBC. Fortunately, we possess an 
efficient tool to substantially reduce the strong finite-size ef- 
fects of PBC by employing screw-periodic boundary condi- 
tions (SBC), as shown recently for the 2D compass model in 
Ref. [9]. Here, we shift the cube L/2 steps in the x-direction 
when leaving the zy-face (plus cyclic permutations), which 
we empirically find to minimize finite-size effects [19]. A nat- 
ural order parameter to detect orbital ordering in the following 
is: 



(3) 



while the complementary quantity D indicates a directional 
ordering of the bond energies: 

D = {l/N)^{E,, - Eyf + (Ey - E,f + {E, - E^f, 

(4) 

which was previously studied in the compass model [7- 
9]. Here, E^^y^^ is the total bond-energy along the 
direction. 

Critical exponents in the Cg model — We start by present- 
ing numerical results for the EgM (1) with SBC by display- 
ing in Fig. 1 the data for the magnetization m and the Binder 
parameter = 1 — (to^)/3(to^)^ as a function of tem- 
perature. Both observables indicate a continuous phase tran- 
sition at about Tc ~ 0.677, in agreement with earlier PBC 
estimates [10, 12]. At Tc we expect Bm{L) to possess only 
coiTections to scaling B„i{L) = + ci^" with uj be- 
ing the correction exponent. We find our best estimate for 
Tc ~ 0.6775(1) and an effective uj ss 1.4 with a large con- 
stant c. 
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FIG. 2. (Color online), eg model: Plot of Xmax and versus L 
in a double logarithmic scale. Estimates for u and rj where obtained 
from a finite-size study using Eq. (5), taking into account corrections 
to scaling. The lines are the corresponding fit curves. 



We now perform a finite-size scaling study to obtain 
the critical exponents. Here, we concentrate primarily on 
the coiTelation length exponent v describing the divergence 
of the correlation length close to the critical point ^ ^ 
|T — Tc| as well as the exponent rj governing the decay 
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of the spin-spin coiTelation function G{r) ^ 
the critical point. We determine these exponents using the 
derivative of the logarithm of the order parameter — 
max{(dlnm/d/3)} [20] and the maximum of the suscepti- 
bility Xmax = max{iV [{w?) — (to)^)} which are known to 
scale with system size L as: [21] 

L^/^l + c^^L-"^), Xmax-i'""(l + CxL-"). (5) 

Using the effective correction exponent lo obtained above 
based on the Binder cumulant, the data fits very well to Eq. (5) 
yielding our estimate v — 0.668(6) for the correlation length 
exponent, see Fig. 2. This value for v would be roughly con- 
sistent with the universality class of the 3D XY universality) 
with vxY — 0.671 [22]. However, an analogous analysis of 
the order parameter correlations at ciiticality - from which we 
obtain rj = 0.15(1) [23] - provides strong evidence for a uni- 
versality class distinct from the 3D XY class, which would 
yield a substantially smaller rjxY ~ 0.038 [6, 22]. Finally, an 
analysis of the exponent a gives a w in agreement with the 
usual hyper-scaling relation. 

Critical exponents in the Cg-clock model — To investigate 
whether the continuous nature of the microscopic degrees of 
freedom T has an impact on the critical properties, we now 
consider a discrete version of Hamiltonian ( 1 ) - one in which 
the vectors T can only point along the six T° ordering direc- 
tions introduced above: 



0/6! 



rii , n 



i-\~e, 



(6) 



Here, E°' {rii , rij ) is the bond energy matrix along the bond 
direction a and n = 0, . . . , 5 denote the six discrete onsite 
states. The similarity of our model to the 6-state (Zg) clock 
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FIG. 3. (Color online), eg-clock model: Orbital order parameter 
m{T) (upper panel) and directional order parameter D{T) (lower 
panel) for different linear system sizes L. Note that both order pa- 
rameters become finite below a common Tc (indicated by the vertical 
line). 

model Hz, = -JE(»,j)T°. • T°^. [24], suggests to tern 
the eg-clock model (EgCLM). Its discrete nature allows 
to study larger systems of up to L = 128. In addition, we 
analyze the directional order parameter D as introduced in 
Eq. (4). In an orbitally ordered state characterized by a finite 
m, D is also finite, however the converse is not true. An il- 
lustrative example is given by the 2D compass model, where 
a gauge-like freedom forbids orbital ordering altogether [25], 
while D orders at finite temperature [7-9]. 

In Fig. 3 we present data for m(T) (upper panel) and D{T) 
(lower panel) for different system sizes. Both m and D ap- 
pear to set in at about the same temperature. In order to con- 
firm the simultaneous onset we have determined the respec- 
tive Binder parameters Bm and Bd (not shown), indicating 
that both transitions take place at a unique critical tempera- 
ture Tc ~ 0.67505(3). This result rules out a scenario of 
a directionally ordered, orbital-disordered intermediate phase, 
and establishes a single transition from a high temperature dis- 
ordered phase to a low temperature orbitally ordered phase. 

Having demonstrated the simultaneity of the two ordering 
phenomena, we now perform a systematic study of the criti- 
cal exponents in the EgCLM. Instead of fitting to Eq. (5), we 
study the finite-size behavior of (running) critical exponents 
obtained on system sizes L and 2L via the relations 

VL=H2)l\n{ml{2L)/m[^{L)), (7) 
77L = 2 - ln(xmax(2L)/x,„ax(^)) /ln(2). (8) 

This allows to visualize finite-size effects directly and should 
give the true exponents for L oo. In Fig. 4 we present 
results for i/^ (upper panel) and rj^ (lower panel). In both 
quantities strong finite-size corrections are evident for the 
EgM and EgCLM, but our results convincingly show that 
different boundary conditions (PBC/SBC) and both the 
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FIG. 4. (Color online). Finite-size scaling of the running exponents 
and r]L calculated from Eqs. (7) and (8) for several models and 
boundary conditions, see legend. Both the eg and the eg-clock model 
display the same critical behavior, which is different from the 3D XY 
universality class (indicated by the dashed lines and similar data for 
the Zs-clock model known to approach 3D XY universality [13]). 

EgM and the EgCLM converge to a single set of exponents: 
1/ w 0.66 and rj w 0.15. These exponents - especially 77 - 
are at variance with the corresponding values of the 3D XY 
universality class. For comparison, we include data for the 
Zg-clock model in Fig. 4, which quickly converges to the 
3D XY exponents expected for this model [13]. Note that a 
similar analysis based on the order parameter D instead of 
m leads to the same v exponent, while the corresponding 
?7d exponent is much larger (w 1.4). This simply follows 
from the assumption that D has no intrinsic critical behavior, 
because then D is driven by m: D ^ m?, resulting in an 
apparently different -q value. 

Emergent U(l) symmetry — In order to shed light on the 
possible emergence of a U(l) symmetry at the critical point 
and the associated behavior of the crossover length scale A 
for T < Tc (as discussed in the context of -perturbed 
XY models [13-16]), we determine the 6-fold anisotropy mg 
of the orbital order m, based on order parameter histograms 
P{r,0) [16]: 

n 1 p2TT 

7716=/ dr ddr^ P{r, 6) cos{6d). (9) 
Jo Jo 

An analysis for the EgCLM analogous to Ref. [16] yields a 
scaling of the crossover length A with the correlation length ^ 
as A C"', with w 1.3 [c.f. Fig. 5(a)]. In the case of a Zg- 
perturbed 3D XY model we find a^^ « 2.2 (c.f. Fig. 5(b), 
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FIG. 5. (Color online). Collapse analysis of me [see Eq. (9)] for the 
eg-clock model (left) and the Ze-clock model (right), based on the 
scaling assumption me ~ L'^^''g{tL^^'-'''^''^) (see Ref. [16]). Best 
collapse parameters ae are indicated in the plot and differ clearly for 
the two models. 



vation, and to understand in more detail the peculiar effects of 
the coupling of real space and order parameter space [6, 27], 
which are at work in these models. Given the broad range of 
systems where models similar to the ones studied here could 
arise (orbital systems in solids [1, 28], Josephson junction 
arrays [5], and artificially engineered systems in optical lat- 
tices [29]), we are optimistic that the peculiar critical prop- 
erties uncovered in the present work can be further explored 
experimentally. 

We thank M. Hasenbusch, G. Misguich, R. Moessner, 
M. Oshikawa, and S. Trebst for useful discussions. The simu- 
lations have been performed on the PKS-AIMS cluster at the 
MPG RZ Garching and on the CalUsto cluster at EPF Lau- 



compatible with Ref. [16]), almost a factor two larger than the 
value we obtain for the EgCLM. 

Compass (t2g) model — Finally we report our results for 
the second orbital-only model of interest here, the t2g model 
in three dimensions [1], defined as: 

= -./ET?Tr+e„, (10) 

where now the degree of freedom T = {T^' ^T^ ^T^-) is a unit 
vector on the sphere S"^, and otherwise the notation follows 
Eq. (1). This model is also called 3D compass model and is a 
straightforward generalization of the 2D compass model stud- 
ied e.g. in [7-9]. An important difference of the model 
compared to the eg model is that orbital order is ruled out 
due to the presence of gauge-like symmetries [25]. Therefore, 
the order parameter D [Eq. (4)] can exhibit a phase transi- 
tion in the absence of orbital ordering. We have simulated 
the full classical t2g model using the same simulation tech- 
nology as for the model, revealing a first order transition 
at Tc « 0.098 from a high-temperature disordered to a low- 
temperature lattice symmetry broken phase indicated by a fi- 
nite value of D. 

Recently the quantum t2g model has been studied using se- 
ries expansions [26], and the absence of a phase transition at 
finite temperature was conjectured. Our findings for the clas- 
sical t2g model provide an alternative explanation as to why 
no (second order) finite-temperature transition was detected: 
due to the first order nature, the transition is intrinsically diffi- 
cult to detect based on series expansions. A detailed analysis 
of the properties of the t2g model will be presented in a forth- 
coming publication [19]. 

Conclusions — We have provided a detailed analysis of 
the critical properties of the finite-temperature ordering tran- 
sitions in eg and t2g orbital-only models. While the <2g model 
exhibits a first order transition, the critical properties of the 
Cg model point towards a distinct universality class, different 
from the standard classes we have encountered so far Further 
theoretical work will be required to shed light on this obser- 
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